Quantum many-body tunneling of attractive Bose-Einstein condensate through 

double asymmetric barrier 
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We study the stability of attractive atomic Bose-Einstein condensate in the anharmonic trap using 
a correlated many-body method. The anharmonic parameter (A) is slowly tuned from harmonic to 
weak and then to strong anharmonicity. For each choice of A the many-body equation is solved 
adiabatically. The use of the van der Waals interaction gives realistic picture which substantially 
differs from the mean-field results. For weak anharmonicity, we observe that the attractive conden- 
sate gains stability with larger number of bosons compared to that in the pure harmonic trap. The 
transition from resonances to bound states with weak anharmonicity also differs significantly from 
the earlier study of Moiseyev etal.[3. Phys. B: At. Mol. Opt. Phys. 37, L193 (2004)]. We also 
study the tunneling of the metastable condensate very close to the critical number N cr of collapse. 
For intermediate anharmonicity, we observe dual tunneling through the two adjacent barriers. We 
also calculate the critical value of A where the left-sided transmission coefficient T; e /t and the right- 
sided transmission coefficient T r i g ht become equal. This is very special feature which has not been 
observed in earlier calculation. For strong anharmonicity we see sharp decrease in the stability of 
the condensate. We observe two separate branches in the stability diagram where we plot N cr with 
tuned A. 

PACS numbers: 03.75.Hh, 31.15.Xj, 03.65.Ge, 03.75.Nt. 



I. INTRODUCTION 

The decay and tunneling of the metastable states is 
an old quantum mechanical problem. However in the 
context of Bose-Einstein condensate (BEC) [_],__[, it re- 
vives new interest in the study of decay and tunneling of 
interacting trapped bosons through a potential barrier. 
Optical traps of finite width can support both the bound 
and resonance states 041] and the macroscopic quantum 
tunneling of a BEC through such finite size barrier is di- 
rectly experimentally observed [1, ■ Theoretically the 
treatment of transport within the mean-field has been 
attempted by many groups ■ Many interesting 

phenomena are observed due to the nonlinear interaction 
term in the Gross-Pitaeveskii (GP) equation. Moiseyev 
et.al. 3 have also studied the transition from resonance 
to bound states of trapped attractive BEC when the in- 
teratomic attraction is increased. Carr et.al. Q have 
studied the macroscopic quantum tunneling in a finite 
potential well in one, two and three dimensions. Time 
dependent GP equation has been used to study the decay 
process of the BEC E jlil. Q__| • In another attempt by the 
group of Heidelberg [ 1 2| , the decay and tunneling dynam- 
ics of few interacting bosons through one-dimensional 
barrier is studied from the first principle and compared 
with the mean-field results. But the quantum many-body 
calculation with realistic interaction for 3D shallow trap 
and the study of stability of attractive BEC with large 
number of bosons is still an open problem. 

In the present work we consider quantum many-body 
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tunneling of the attractive BEC near the criticality and 
trapped in a quadratic with a quartic confinement. The 
quartic term takes care of the shallow Gaussian poten- 
tial of finite width. The external potential is modeled as 
V(r) = ^mu> 2 r 2 — Ar 4 . This represents the optical trap 
used in many experimental BEC [HI, [l4[ . This type of 
study is important for various reasons. First, near the 
critical point, as the condensate becomes highly corre- 
lated the mean-field GP equation does not produce good 
results. This is pointed earlier that the choice of contact 
interaction in 3D attractive condensate does not repre- 
sent the true interatomic interaction [l6[ . So near the 
critical point where the BEC has high probability of tun- 
neling through the adjacent barrier one needs to apply 
an approximate many-body technique which takes care 
the interatomic correlation. Second, the choice of van der 
Waals interaction in the many-body theory with a short 
range hard sphere and a — ^f- tail gives the realistic pic- 
ture. In the many-body effective potential the metastable 
BEC is now bounded by double Gaussian-like barriers of 
different height in two sides. Thus the BEC suffers tun- 
neling through both the barriers simultaneously. This is 
distinctly different from the tunneling in the mean-field 
approach where the BEC suffers tunneling through the 
right side barrier which is the effect of finite trap size 
only. Although a preliminary attempt has been made 
in our earlier published work in this direction [I?} but 
it needs further study for the following reasons. In our 
earlier study we have considered tunneling for a fixed A 
parameter. In the present work we continuously tune A 
and study the tunneling process adiabatically. This needs 
to control the height of the confining potential well in a 
controlled fashion. As the trapping potential is imposed 
optically by laser, in the routine experiments of BEC, the 
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height of the external potential is controlled by reducing 
the laser intensity. Thus our present study is quite im- 
portant as it could be obtained in the laboratory with 
present day set-up. 

In our present calculation we choose 85 Rb condensate 
in the anharmonic trap. We study the various tunnel- 
ing phenomena for different choices of anharmonic pa- 
rameter A: weak, intermediate and strong. We observe 
new physics which substantially differs from mean-field 
picture. For very weak anharmonicity the stability of 
the attractive Bose condensate is substantially improved 
whereas for very strong anharmonicity the condensate 
collapses quickly. For intermediate anharmonicity we ob- 
serve dual tunneling of the condensate through the adja- 
cent barriers. 

The paper is organized as follows. In Section II, we 
introduce the many-body calculation with the correlated 
harmonic potential basis. Section III discuses the numer- 
ical results and Section IV concludes the summary of our 
work. 



II. METHODOLOGY 

A. Many-body calculation with correlated 
potential harmonic basis 

The earlier theoretical studies on attractive BEC in the 
harmonic trap used the mean-field approximation which 
results the Gross-Pitaeveskii (GP) equation for contact 
5- interaction [l8| . As the total condensate wave function 
is taken as the product of single particle wave functions, 
the effect of interatomic correlation is completely ignored. 
However specially for the attractive condensate, as the 
atoms come closer and closer, the central density becomes 
high and the condensate becomes highly correlated near 
the critical point. Naturally the interatomic correlation 
can no longer be ignored and one needs a full quantum 
many-body calculation which takes care of the effect of 
interatomic correlation. 

In our present study we solve the many-body 
Schrodinger equation by potential harmonic expansion 
method (PHEM), which basically uses a truncated two- 
body basis set which keeps all possible two-body correla- 
tion [l9[ and we go beyond the mean-field approximation. 
The potential harmonic expansion method with an addi- 
tional short range correlation function, called CPHEM, 
has already been established as a very useful technique 
for the study of attractive BEC 0, lH H3 . Here we 
describe the methodology briefly for the interested read- 
ers. Details are found in our earlier work [22|, [23|. The 
Hamiltonian for a system of A = (N + 1) identical 
bosons (each of mass m) interacting via two-body po- 
tential V(fij) = V(fi — fj) and confined in an external 
trap (which is modeled as a harmonic potential with a 



quartic term) has the form 

H = E v » + E W?o + E nn - r,y (i) 

i— 1 i— 1 i,3>i 

After elimination of the center of mass motion and using 
standard Jacobi coordinates [HI, HH, [25| , the Hamiltonian 
describing the relative motion of the atoms is given by 




i=i 



Vi n t is the sum of all pair-wise interactions expressed in 
terms of the Jacobi vectors. It is to be noted that Hyper- 
spherical harmonic expansion method (HHEM) is an ab- 
initio tool to solve the many-body Schrodinger equation 
where the total wave function is exp anded in the com- 
plete set of hyperspherical basis [2J|. Although HHEM 
is a complete many-body approach which includes all 
correlations, due to large degeneracy of the HH basis, 
HHEM can not be applied to a typical BEC which con- 
tains few thousands to few millions of atoms. However in 
the context of experimentally achieved BEC, as the in- 
tcrparticle separation is very large compared to the range 
of interatomic interaction, we can safely ignore the effect 
of three-body and higher-body correlation and can keep 
only the two-body correlation. This is perfectly justified 
for dilute BEC where the probabilities of three and higher 
body collision is negligible. It permits us to decompose 
the total wave function "J into two-body Faddccv com- 
ponent for the interacting (ij) pair as 

A 

*= E <Mw) • (3) 

i,j>i 

It is worth to note that 4>ij is a function of two-body sepa- 
ration (fij ) only and also includes the global hyperradius 

r, which is given 

by r = v/E£LiC?- Thus thc effect of 
two-body correlation comes through the two-body inter- 
action in the expansion basis. <f>ij is symmetric under P^- 
for bosonic atoms and satisfy the Faddeev equation 

A 

[T + V trap - E R ] fa = -V(rij ) <S>u (4) 

kl>k 

where T is the total kinetic energy. Operating . >{ 
on both sides of equation (4), we get back the original 
Schrodinger equation. In this approach, we assume that 
when (ij) pair interacts, the rest of the bosons are in- 
ert spectators. Thus the total hyperangular momentum 
quantum number as also the orbital angular momentum 
of the whole system is contributed by the interacting pair 
only. Next the (ij')th Faddccv component is expanded in 
the set of potential harmonics (PH) (which is a subset 
of HH basis and sufficient for the expansion of V (?%•)) 
appropriate for the (ij) partition as 

M?ii>r)=r-^ ) ^V&«{il%)u l K (r) ■ (5) 

K 
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f2)J denotes the full set of hypcrangles in the SAT- 
dimensional space corresponding to the (ij)th interacting 

pair and P^k+A^n) i s called the PH basis. It has an an- 
alytic expression: 



1 2K+iy 1L N 



i) {N) P l 2 K+l ^)y (D-i)- 



D = 3 AT, 
(6) 

3^o {D — 3) is the HH of order zero in the (3 AT — 3) dimen- 
sional space spanned by ...,£jv-i} Jacobi vectors; <j) 
is the hypcranglc given by = r coscj). For the remain- 
ing (AT — 1) noninteracting bosons we define hyperradius 
as 



Pij 



N-l 



r suit 



(7) 

such that r 2 = rfj + pf^ and r represents the global hyper- 
radius of the condensate. The set of (3A~ — 1) quantum 
numbers of HH is now reduced to only 3 as for the (AT — 1) 
non-interacting pair 

h = h = - = In-i = 0, (8) 
mi = ni2 = ... = ra/v-i = 0, (9) 
n,2 = n 3 = ...njv-i = 0, (10) 

and for the interacting pair In = I, = m and 

tin = K. Thus the 3AT dimensional Schrodinger equation 
reduces effectively to a four dimensional equation with 
the relevant set of quantum numbers: hyperradius r, or- 
bital angular momentum quantum number I, azimuthal 
quantum number m and grand orbital quantum number 
2K + I for any AT. Substituting in Eq(4) and projecting 
on a particular PH, a set of coupled differential equation 
(CDE) for the partial wave u l K {r) is obtained 



77i dr' 2 



W) + ^{C{£ + 1) 
+4K(K + a + p + 1)} - E R ] U m {r) 



(11) 



+ YjKiV K K>(r)fK>iU K >i(r) =0 



K 



U K i = fKiu l K (r), 



Q 



where C = I + 24j=5, 
= 1 + 1/2. 

fxi is a constant and represents the overlap of the PH 
for interacting partition with the sum of PHs correspond- 
ing to all partitions (25[. The potential matrix element 
Vkk > (r) is given by 



V KK '(r) = 



(12) 



Introduction of additional short range 
correlation 



As pointed earlier in the mean-field GP equation the 
two-body interaction is represented by a single parame- 
ter, the s-wave scattering length a sc only. It disregards 



the detailed structure. The presence of essential singu- 
larity as r — > for the attractive contact <5-intcraction 
makes the Hamiltonian unbound from below. So in our 
present many-body calculation, we use a realistic inter- 
atomic potential like the van der Waals potential with 
an attractive — tail at large separation and a strong 
short range repulsion. The inclusion of detailed structure 
in the two-body potential with the short range repulsive 
core needs to include an additional short range correla- 
tion in the PH basis. This short range behavior is repre- 
sented by a hard core of radius r c and we calculate the 
two-body wave function rj(rij) by solving the zero-energy 
two-body Schrodinger equation 



h 2 1 d 



di](rj 
dru 



+ V(r ij )r ) (n j )=0- (13) 



This zero-energy two-body wave function ^(ry) is a good 
representation of the short range behavior of 4>ij as in the 
experimental BEC, the energy of the interacting pair is 
negligible compared with the depth of the interatomic po- 
tential. It is taken as the two-body correlation function 
in the PH expansion basis. The value of r c is obtained 
from the requirement that the calculated a sc has the ex- 
pected value [i^l . We introduce this as a short-range 
correlation function in the expansion basis. This also 
improves largely the rate of convergence of the PH basis 
and we call it as correlated Potential Harmonic expansion 
method (CPHEM). We replace Eq(5) by 



E^ra^w^)- ( 14 ) 



K 



and the correlated PH (CPH) basis is given by 



[V. 



+1 (^iV ^ > 



correlated 



2K+1 



(rfipHm), (is) 



3A-8 anc j The correlated potential matrix Vkk 1 ( r ) is now given by 




Vkk\t) 

/+ 1 {P^(z)^(r v /l^)P^(z) 7? (r v /l^)^(z)}dz 

(16) 

Here P K ^(z) is the Jacobi polynomial, and its norm and 
weight function are h^f and Wi(z) respectively [27^ 

One may note that the inclusion of rj(rij ) makes the PH 
basis non-orthogonal. One may surely use the standard 
procedure for handling non-orthogonal basis. However in 
the present calculation we have checked that Tjfcj) dif- 
fers from a constant value only by small amount and 

the overlap (^K+i^^^l^x+iC^iV^) 7 ?^)) is q uitc 
small. Thus we get back the Eq(ll) approximately when 
the correlated potential matrix is calculated by Eq(16). 
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III. RESULTS 

A. Choice of two-body potential and calculation of 
many-body effective potential 

The interatomic potential has been chosen as the van 
der Waals potential with a hard core of radius r c , viz., 



for 



< 



and 



c 6 



for 



> r c , C 6 is 



known for a specific atom and in the limit of Cq — > 0, 
the potential becomes a hard sphere and the cutoff ra- 
dius exactly coincides with the s-wave scattering length 
a sc . For our present study we consider 85 Rb atoms in the 
JILA trap where the controlled collapse experiment has 
been observed [28|, [29[. By utilizing the Feshbach res- 
onance one can effectively tune the negative scattering 
length and can study the stability in controlled fashion. 
In our choice of two-body potential we tune r c to repro- 
duce the experimental scattering length. As we decrease 
r c , a sc decreases and at a particular critical value of r c 
it passes through — oo to oo [Hj]. For our present calcu- 
lation we choose 85 Rb atoms with C* 6 = 6.4898 x 1CT 4 
o.u. [l[ and a sc ~ —1.832 x 10~ 4 o.u., which is one of the 
choices of scattering length a sf in the controlled collapse 
experiment of Roberts et.al. [H, [29|. With this Vfaj) 
we solve the zero-energy two-body Schrodingcr equation 
and tune r c to obtain correctly a sc = —1.832 x 10~ 4 o.u. 
We choose r c such that it corresponds to the zero node in 
the two-body wave function. The chosen parameter for 
our calculation is r c = 1.3955 x 10~ 3 o.u.. With these set 
of parameters we solve the coupled differential equation 
by hypersphcrical adiabatic approximation (HAA) f30j ] . 
In HAA, we assume that the hyperradial motion is slow 
compared to the hyperangular motion and the poten- 
tial matrix together with the hypercentrifugal repulsion 
is diagonalized for a fixed value of r. Thus the effec- 
tive potential for the hyperradial motion is obtained as 
a parametric function of r. We choose the lowest eigen 
potential ujq (r) as the effective potential in which the con- 
densate moves collectively. The energy and wave function 
of the condensate are finally obtained by solving the adi- 
abatically separated hyperradial equation in the extreme 
adiabatic approximation (EAA) 



— -7-0 + ^w) - e r 



Co(r)=0, 



(17) 



subject to approximate boundary condition on (a(r). For 
our numerical calculation we fix I — and truncate the 



CPH basis to a maximum value K = K„ 



requiring 



proper convergence. In Fig.l we plot the many-body ef- 
fective potential u>o(r) as a function of hyperradius r for 
N = 1000 atoms in the pure harmonic trap (A = 0). For 
N < iV cr (~ 2385), the condensate is metastable and is 
associated with a deep and narrow attractive well (NAW) 
on the left side. For r — > 0, there is a strong repulsive 
wall which is the reflection of the hard core van der Waals 
interaction. For N less than the critical number N cr a 
metastable region (MSR) appears for larger r. An inter- 




1.5 



(a) 




(h) 



FIG. 1: (Color online) Plot of effective potential u)o(r) of an 
attractive BEC with N = 1000 in a pure harmonic trap. In 
panel (a) the deep narrow attractive well NAW is ploted while 
in panel (b) the metastable region MSR is presented. 



mediate barrier (IB) separates the NAW from the MSR. 
The panel (a) of Fig. 1, the NAW together with the repul- 
sive core is shown. The IB and MSR have been shown in 
panel (b) of Fig. 1. With increase in N, we observe that 
the height of IB decreases, together with a decrease in 
the difference (Aw) between the maximum of IB (uj max ) 
and minimum of MSR (u) m in) and the NAW starts to be 
more negative and narrower. As N — > N cr , uj max and 
u m in coincides, the MSR disappears which corresponds 
to the collapse of metastable BEC. For N > N cr , there 
will be only the NAW and the condensate will settele 
in the deep attractive well. In the attractive (5-function 
potential, as r — > 0, the effective potential rapidly goes 
to —00 which is called a pathological singularity. Here 
the many-body picture is in sharp contrast with the GP 
mean-field picture with contact interaction. 
The advantages of the CPHEM are as follows. 

i) As for the attractive BEC, the condensate becomes 
highly correlated, the choice of CPH basis (which 
keeps all possible two-body correlations) is per- 
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fectly justified for the description of correlated di- 
lute BEC. 



ii) By using the HAA we basically reduced the 
3A dimensional problem into an effective one- 
dimensional problem in hyperradial space which 
provides both qualitative and quantitative picture 
of the system. As the effective quantum numbers 
are always four (for any N) we can run our many- 
body code for vary high values of TV. 

iii) The many-body effective potential strongly differs 
from the mean-field potential. In GP, the choice 
of contact (^-interaction in the two-body potential 
gives rise to the pathological singularity in the ef- 
fective potential. Thus the study of post-collapse 
scenario of the attractive condensate is beyond the 
scope of mean-field theory. Whereas the presence 
of short-range hard core in the van der Waals inter- 
action does not only remove the singularity, it gives 
the realistic scenario and can describe the forma- 
tion of atomic cluster after collapse. 



B. Weak anharmonicity (A < 9.36 x 1CT 6 o.u) 

In our earlier work we have succesfully applied 
CPHEM to describe controlled collapse experiment and 
calculated the stability factor more accurately than the 
mean- field results [l6T ] . Thus the main aim of our present 
work is to study the stability of 85 Rb atoms with fixed 
interaction when the shallowness of external trapping po- 
tential (A) is continuosly tuned. Experimentally the tun- 
ing is done by controlling the intensity of laser beam, 
which provides the optical trap. The many-body ef- 
fective potential for weak anharmonicity is plotted in 
Fig. 2(a) for N = 2395 which is very close to critical- 
ity N cr = 2400. We slowly change A, from very close to 
harmonic to A = 9.35 x 10~ 6 o.u. and numerically solve 
the many-body equation adiabatically for each choice of 
A. For such weak anharmonicity, the mctastablc conden- 
sate in the MSR region is tightly bounded by the high 
intermediate barrier on the right side (RIB) and does 
not suffer any tunneling through RIB. Whereas the inter- 
mediate barrier on the left side (LIB) is of finite height 
and the condensate suffers tunneling through the LIB 
only. The LIB is again associated with a deep negative 
well (NAW) on the left side as shown in the Fig. 1 and 
the condensate after tunneling settles down in the NAW 
and forms cluster. By increasing A gradually, we observe 
the greater stability of the attractive condensate as the 
metastable region (MSR) becomes deeper. The stability 
of the condensate with gradually increasing values of A 
is also clear from Fig. 2(b) where we plot the wavefunc- 
tion for various choices of A. The wave function expands 
with increasing A. From the condensate wave function 



we calculate the condensate radius r av as 



G 



1/2 < r 2 > l/2 
V2A ' 



(18) 



where R is the center of mass coordinate. And we plot 
it in Fig. 2(c). We observe that r av increases gradually 
with increase in A. We also calculate the WKB tunneling 
probability T as 



T = cxp(-2 / V 2 M0 -E] dr) 
J n 



(19) 



where the limits of integration r± and r2 are the inner 
and outer turning points of the intermediate barrier on 
the left (LIB) of uio(r), E is the energy of the metastable 
condensate. We plot it as a function of A in the Fig. 3. It 
nicely presents that for larger A the condensate gets bet- 
ter stability and less tunneling. Thus tuning A gradually 
from very close to harmonic to weak anharmonicity we 
observe increasing stability of the condensate. Thus sim- 
ply changing shape of the external trap one can produce 
highly stable condensate even for the attractive bosons. 
Whereas with decreasing A gradually as the transmission 
coefficient T increases, the system will collapse for a cer- 
tain number of atoms. This is similar to the collapse of 
attractive BEC in pure harmonic trap. In our present cal- 
culation with the choice of parameters we observe greater 
stability till A < 9.36 x 10~ 6 o.u.. 



C. Intermediate anharmonicity 

(9.36 x 1(T 6 < A < 9.53 x 1(T 6 ) 

This is very interesting situation when the condensate 
in the MSR is bounded by two barriers (LIB and RIB) 
of finite width and may suffer tunneling from the MSR 
through both the barriers simultaneously For better un- 
derstanding we plot the many-body effective potential in 
Fig. 4 for some specific choices of A within the above 
range and for fixed A=2665. For A = 9.36 x 10~ 6 , 
the height of LIB is smaller than the height of RIB. 
Thus although the condensate will suffer dual tunneling, 
however the transmission probability through the LIB is 
heigher compared to that through RIB. Increasing A to 
9.40 x 10 -6 , we observe that the height of LIB and RIB 
are comparable and the tunneling probability of the con- 
densate through LIB and RIB are of equal weightagc. 
Increasing A further, we observe that although the LIB 
remains of significant height, RIB drastically falls. Natu- 
rally the condensate will suffer larger tunneling through 
the RIB. As stated before there is a large deep negative 
well on the left side of the LIB, which is not shown in 
the figure for better presentation. Thus when the con- 
densate leaks through the LIB, it is not lost. The highly 
correlated bosons will form cluster within the deep well. 
Whereas after tunneling through the RIB the interact- 
ing bosons will be traplcss and will behave as uniform 
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Bose gas. The calculated WKB tunneling probability T 
is plotted as a function of A in this intermediate region 
in Fig. 5. The double branches definitely signify the 
dual possibility of tunneling through the double barrier. 
However unlike the double Gaussian barrier, we have two 
asymmetric barrier of different height. In Fig. 5, we ob- 
serve that the tunneling through LIB (T; e / t ) gradually 
decreases with A which is similar to the case of weak 
anharmonicity as discussed previously. Whereas the 
transmission probability through RIB (T rig i lt ) sharply in- 
creases with A which is similar to the case of strong an- 
harmonicity (described in the next section). We see that 
Ti e f t and T r i g ht cross at some value of A. The crossing 
point is shown in the inset of Fig. 5. The correspond- 
ing value of A is 9.398 x 10~ 6 o.u. At this critical A 
(^critical), Ti e f t becomes equal to T right , which signifies 
that the two asymmetric barriers become two symmet- 
ric Gaussian barriers and the condensate will suffer equal 
tunneling. Thus our observation nicely demonstrate that 
by simply controlling the external laser beam, one can 
control the nature of anharmonicity and stability of the 
condensate. The value of X cr iticai depends on the scatter- 
ing length a sc . Thus taking different choices of a sc from 
controlled collapse experiments of Roberts [28| , one can 
calculate different critical values of A. 



D. Strong anharmonicity (A > 9.53 x 10 6 ) 

This is the other end of anharmonicity where the con- 
densate is tightly bound in the MSR by a high barrier 
on the left side (LIB) and suffers no tunneling through 
it. But the RIB is of finite height and the condensate 
suffers tunneling through RIB only. With increase in A 
value, the barrier height of RIB quickly decreases which 
enhances the tunneling probability rapidly. The effective 
potential is plotted in the Fig. 6(a). Unlike the previous 
case of weak anharmonicity, here the MSR has no deep 
well on the right side. So after tunneling the condensate 
will behave as trapless, uniform Bose gas. The feature 
is also clear from Fig. 6(b) where we observe that the 
condensate wave function is associated with oscillation 
on the right side. For better clarity the oscillatory part 
of the wave function is presented in the inset of Fig. 6(b) . 
It clearly shows that the condensate is basically confined 
in the MSR but a part likages through the adjacent RIB. 
With increasing A the oscillation also increases. The cor- 
responding WKB tunneling probability is presented in 
Fig. 7 which shows sharp increase in T even with small 
change in A which indicates the quicker collapse of the 
condensate in a strong anharmonic trap. 



E. Collapse of atomic BEC in anharmonic trap 

Here we study the collapse of attractive Bosc-Einstein 
condensate in the shallow trap which is produced with 
an additional quartic anharmonic potential with a har- 



monic potential. In the BEC experiments of attractive 
condensate [H, HI] , it is observed that the condensate is 
metastablc for a finite number of atoms. This is called 
critical number (N cr ) and the condensate collapses when 
N > N cr . In the controlled collapse experiment of 85 Rb 
condensate, the N cr is calculated both by the GP mean- 
field theory [l8j and ab-initio quantum many-body calcu- 
lation [16| . The critical stability factor is more accurately 
calculated in the correlated many-body approach [26| and 
is closer to the experimental value compared to the mean- 
field result. Thus in the present study we are interested 
to calculate N cr when the external trap is slowly tuned 
from harmonic to weak anharmonic and then strong an- 
harmonic. In Fig. 8 we plot the many-body effective 
potential for both harmonic and with weak anharmonic- 
ity. For A = 0, the trap is purely isotropic, harmonic 
and the condensate collapses when N > N cr . In Fig. 8 
we observe that for A = 0, the intermediate barrier on 
the left side (IB) just vanishes for N = 2385 and the 
metastable condensate confined in the MSR settles down 
in the negative deep well on the left side of the intermedi- 
ate barrier (IB) (not shown in the figure). It corresponds 
to collapse and the corresponding critical number for col- 
lapse is N cr = 2385. Now with the harmonic potential if 
we add anharmonic potential with small anharmonicity 
(A = 1.0 x 10~ 6 ), the MSR starts to develop and the inter- 
mediate barrier reappears. For this A, the condensate col- 
lapses for N cr = 2400. Thus by increasing A by very small 
amount, the critical number increases by 15 atoms and 
gives better stability of the condensate. It corresponds 
to the metastable condensate which is confined in MSR. 
Thus the stability is enhanced by controlling A. The 
large r behavior of the eigen potential almost remains 
unchanged. Thus by increasing A very slowly we observe 
the condensate becomes more stable and N cr gradually 
increases with A. It is shown in Fig. 9, which nicely 
describes that N cr smoothly increases with A, when A 
is small. In the earlier study of Moiseyev et.al. [|[ the 
transition from resonance to bound state has been stud- 
ied when the attractive interaction gradually increases. 
However to study such transition, an additional negative 
offset (Vo) was required to prevent collapse. For Vq = 0, 
the authors found in the ref. @ that resonance never 
turns into bound state as the nonlinearity becomes more 
negative, the condensate collapses. However due to the 
use of hard core van der Waals potential we do not re- 
quire any such external offset to create the metastable 
condensate. Thus the increasing stability and the transi- 
tion from resonance to quasi-bound states for very weak 
anharmonicity strongly differs from the earlier observa- 
tion of Moiseyev et.al. . However with further increase in 
A, we enter in the zone of intermediate anharmonicity and 
strong anharmonicity which shows that N cr now sharply 
decreases with A. Thus we observe two branches in Fig. 
9. Wc very carefully studied and observed that these 
two branches never meet. The magnified graph near the 
point of bifurcation is also presented in the inset of Fig. 
9. 
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IV. CONCLUSION 

In our present work, we have studied the decay and 
tunneling of attractive Bose Einstein condensate in the 
finite trap by approximate but ab-initio many-body cal- 
culation. Main attention has been paid in the study of 
stability of attractive condensate when the effective trap 
height is tuned from very close to harmonic to weak an- 
harmonicity and then to strong anharmonicity. The use 
of correlated PH basis and the realistic van der Waals 
interaction correctly describes the many-body tunnel- 
ing process and gives the realistic picture. Our result 
strongly differs from the earlier study of mean-field re- 
sults through onc-dimcnsional barrier. The greater sta- 
bility of attractive Bose gas in the anharmonic trap with 
very weak anharmonicity and the dual tunneling of the 
attractive condensate for intermediate anharmonicity are 
of special interest. The use of hard core potential in the 
two-body interaction can always describe the post col- 
lapse scenario of the attractive condensate in the many- 
body picture. No additional negative offset potential is 
required at r = to observe collapsc-to-mctastable tran- 
sition. Such attractive condensates in the finite trap are 
created in the laboratory when the height of the confining 
potential is controlled by controlling the laser intensity in 
the optical trap or by controlling the electric current in 
magnetic trap. Thus the attractive BEC can be created 
in a shallow trap of chosen height and once it is material- 



ized several experiments can be done. Thus our present 
many-body calculation is quite relevant experimentally 
and calculated observables like the tunneling probability 
T and critical value of A parameter, N cr can be easily 
observed in the present day experiments. The choice of 
scattering length a sc is also taken from the experiment 
of 85 Rb atoms in JILA. 

This is the first such many-body calculation where we 
compute resonances and decay of many-boson systems 
in three-dimension beyond the mean-field. A natural ex- 
tension of the present work is to study the dynamics and 
to calculate the non-escape probability and various co- 
herence properties for many-boson systems in 3D, which 
requires the solution of full time dependent many-body 
equation. 
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FIG. 2: (color online) a) Plot of effective potential uJo(r) of 
the condensate with N = 2395 for different A; b) Plot of cor- 
responding condensate wave functions and c) Plot of the aver- 
age condensate S1Z6 Tav for TV = 2395 against anharmonicity 
parameter A. 
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FIG. 3: (color online) Plot of Transmission coefficient T vs 
the anharmonicity parameter A for N = 2395. 




FIG. 4: (color online) Plot of effective potential u>o(r) for 
different anharmonicity parameter A for N = 2665. 




FIG. 5: (color online) Plot of Transmission coefficient T vs 
anharmonicity parameter A. The inset shows the crossing 
point of Tle/t and T r i g ht- 
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FIG. 6: (color online) a) Plot of effective potential uJo(r) of 
the condensate with N = 160 for different A; b) Plot of corre- 
sponding wave functions of the condensate . In the inset the 
oscillatory part of the condensate wave function is shown. 
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FIG. 7: (color online) Plot of Transmission coefficient T vs the 
anharmonicity parameter A for a condensate with N = 160. 
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FIG. 8: (color online) Plot of effective potential u> for different 
weak anharmonicty A for a condensate with N — 2385. 
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FIG. 9: (color online) Plot of critical number N cr vs the an- 
harmonicity parameter A. 



